Exponential velocity tails in a driven inelastic Maxwell model 
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The problem of the steady-state velocity distribution in a driven inelastic Maxwell model of shaken 
granular material is revisited. Numerical solution of the master equation and analytical arguments 
show that the model has bilateral exponential velocity tails {P{v) ~ e"'"'''^^), where D is the 
amplitude of the noise. Previous study of this model predicted Gaussian tails {P{v) e""" ). 



Recently, granular systems have been intensively stud- 
ied ||. One reason of such an interest is the fact that 
in many respects these systems are highly unconventional 
and even exotic As a result, very often even stan- 
dard laws of statistical physics require certain modifica- 
tion when applied to granular systems. 

As an example, let us consider the Maxwell-Boltzmann 
law for velocity distribution P{v) of atoms or molecules in 
a gaseous state, which states that P{v) ~ e""" . Under 
certain experimental conditions, a granular system can 
be considered as a gas and a natural question is what 
is its velocity distribution. Numerous theoretical and ex- 
perimental studies do not provide a simple answer to this 
question. On the contrary, they show that P{v) depends 
on certain details of the system as for example how en- 
ergy is transferred by a thermostat into the system, in 
order to balance the energy loss during inelastic colli- 
sions. Theoretical works have shown that P{v) might be 
of the form e'"''^ g |5|, e-"^""^" § | or e""^" @. Un- 
der certain conditions experiments show clear deviations 
from the Maxwell-Boltzmann law, but it is still rather 
difficult to decide what is the form of P{v) in real gran- 
ular systems 

Very often granular systems are described using the 
so-called inelastic hard-sphere models, which could be 
then analysed using corresponding Boltzmann equations. 
These off-lattice, two- or three-dimensional systems are 
however very difficult to study, especially when we want 
to explore large-velocity regions of the phase space. A 
possible alternative is to construct simplified models for 
which more accurate calculations are possible. One class 
of such models are Maxwell models, for which the colli- 
sion term in the corresponding Boltzmann equation is ve- 
locity independent. Recent studies show that for Maxwell 
models P{v) might take cither exponential (e"'"') jlo| or 
Gaussian (e"-^"") form |ll|, 

A particularly interesting Maxwell model was proposed 
by Ben-Naim and Krapivsky (BN-K) . These authors 
presented an elegant solution of the master equation of 
their model using a Fourier transform method. However, 
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their conclusion that P{v) has a Gaussian decay is based 
on certain approximation whose validity is difficult to 
assess. Namely, they infer the large- velocity behaviour of 
P{v) from the \ow-k behaviour of the Fourier transform 
P{k). Since BN-K model is one of the very few models 
for which exact or numerical but very precise calculations 
can be made, it would be desirable to clarify the validity 
of this approach. 

In the present paper we reexamine the BN-K model. 
Analysing numerically the solution of the master equa- 
tion of the model we obtain the velocity distribution 
P{v). Asymptotically (v — > oo), this quantity shows 
bilateral exponential tails {P{v) ^ e""'"') and such a 
behaviour is seen over more than 10 decades. We also 
present analytical arguments which support existence of 
bilateral exponential tails in this model. From our anal- 
ysis it follows that a ~ where D is the amplitude of 
noise which simulates the input of energy into the system. 

To introduce the model, we consider a collection of 
particles that are characterized by a single parameter, 
their velocity v. In this model positions of particles are 
not specified hence the model neglects any spatial correla- 
tions. Particles undergo two-body inelastic collisions that 
change their velocities according to (wi,W2) — > {^[,'^2), 
where 



7 

1-7 



1-7 



(1) 



and 7 is the inelasticity parameter (0 < 7 < 1). Particles 
which participate in a collision are chosen at random. In 
addition to the collision rule (p, particles are subjected 
to the uncorrelated white noise of strength D. Existence 
of noise ensures that the model has a well-defined non- 
trivial steady state. In the steady-state the velocity dis- 



tribution P{v) satisfies the following equation 12 



-DP"{v) = -P{v) 
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where double prime " denotes the second derivative with 
respect to velocity. To solve eq. BN-K introduced 
the Fourier transform of the velocity distribution P{k) = 
J dve^^"" P{v), which in the steady state satisfies 



(1 -I- Dk^)P{k) = P[7fc]P[(l - 7)fc]. 



(3) 
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This equation admits the following solution 

oo i 

P{k)^\{Y{[l + l'''{l-lf^'-''^De]. (4) 
After some transformations eq. (|^) can be written as 



P{k) — exp 



E 



(5) 



where a„(7) = 1 — (1 — 7)" — 7". To obtain a large 
velocity behaviour of P{v) BN-K truncate the series in 
eq. keeping only the first, quadratic in k, term 



P{k) ~ exp 
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Subsequently, taking the inverse Fourier transform of this 
Gaussian function they obtain P{v) ^ exp(-=^), where 

vl=2D/a2{i). 

In our opinion, the above procedure of truncating the 
series and inverting the resulting function is not well jus- 
tified. The approximation (^) agrees with the exact ex- 
pression (||) only up to the fc^ term of the Taylor ex- 
pansion. Equally consistent approximation can be given 
as 



P{k) 



k^' 



(7) 



where c = y^^- But the inverse Fourier transform 

of (|^) is |exp(— c|u|), which is qualitatively different than 
the Gaussian decay obtained by BN-K. Both (^) and (|^) 
are consistent with the exact solution (||) up to the k^ 
order but differ at the higher orders. Apparently, these 
higher order terms qualitatively affect the large-velocity 
limit of P{v). Consequently, without more detailed ar- 
guments such approximations are not justified. Indeed, 
as we will show below, P{v) has a bilateral exponential 
decay but with a different coefhcient c. 

To check the validity of the BN-K approach we cal- 
culate P{v) numerically as an inverse Fourier transform 
of P{k). One way to compute P{k) is to evaluate the 
infinite sums in the logarithm of eq. (^. However, it 
turns out that using eq. @ leads to a better precision. 
Let us notice, however, that the series in eq. (j^) con- 
verges but only for |fc| < fee = To calculate P{k) for 
|fc| > kc we can then use eq. (y|), provided that jk and 
(1 — "f)k fall within the range of convergence. If not, we 
have to calculate Pijk) and P((l — 7)fc) referring again 
to eq. (H). Implementing this recursive procedure, we 
calculated P{k) and then using the Fast Fourier Trans- 
formation algorithm we obtained P{v). Our results for 
D = 1 are shown in Fig. |l|. This figure clearly shows that 
P{v) decays exponentially. For 7 = 0.5 and 0.8 such a 
behaviour is seen for more than 10 decades. For 7 close 
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FIG. 1: The velocity distribution P{v) as a function of v 
calculated for D = 1. Due to symmetry only v > Q part 
is shown. The thick line represents the exponential function 
exp(— w). 



which asymptotically is replaced by the bilateral expo- 
nential decay (exp(— c|t;|)). Moreover, as shown in Fig. |l| 
the asymptotic slope c is independent of 7 and is ap- 
proximately one. As we will show below this slope only 
depends on the amplitude of noise D. Before doing that 
let us notice that the D-dependence of the velocity dis- 
tribution P(y) can be easily inferred from the fact that 
D and k enter its Fourier transform only through Dk^ 
terms (see eqs. and (^)). From this property one can 
easily obtain that 



(8) 



where we explicitly indicated the dependence on the noise 
D. 

In the following we provide some analytical arguments 
showing that P{v) decays exponentially. Let us assume 
that in eq. (||) the second term (gain) can be neglected. 
Then eq. (||) simplifies to DP"{v) = P{v) and the nor- 
malizable solution reads 



P{v) 



(9) 



to unity there is a exp(- 
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decay for small velocities 



Let us notice that the exponential part is 7-independcnt 
(since 7 enters the master equation only through the ne- 
glected gain term). Moreover, such a solution is in a very 
good agreement with numerical calculations (see Fig. |^). 

For more general models it is known that neglecting the 
gain term is justified for v —^ 00 when the resulting solu- 
tion decays faster than exponentially ^, . Our solution 
(^) is thus a marginal case. However, we can show that 
for our model the gain term in the limit v 00 indeed 
can be neglected. First, let us evaluate the gain term for 
the solution (^. Elementary integration for Z? = 1 and 
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One can see that since < 7 < 1, the gain term decays 
exponentially with v but faster than the solution (^) (i.e., 
a non- neglected loss term in (||)). We expect that eq. (||) 
is only an asymptotic (v ^ 00) solution of the master 
equation (||). Thus, for velocities u ~ 0(1) the distribu- 
tion P{u) deviates from the asymptotic form (^. This 
will modify the integral ( p^ ) but only in the vicinity of 
u — and u — ^. As we argue below, such a modifica- 
tion of P{u) still leads to the gain term decaying faster 
than the loss term. 

Indeed, the contribution around u = is a prod- 
uct of P{u) and of the exponential term exp( ~^~^" ) . 
Thus, for large v a modification of P{u) around u = 
will change our estimation ( p^ ) but only at the order of 
exp(j^^). Similarly, one can show that a modification 
of P{u) around u — will change ( p^ by a factor of 
the order of exp(^). Consequently, the gain term again 
decays faster than the solution P{v), which justifies its 
neglect. Let us also notice that when 7 approaches or 
1 the model becomes energy conserving and the bilat- 
eral exponential distribution is no longer expected to 
hold. But it is easy to notice that in such cases the decay 
of the gain term matches the decay of the solution P{v) 
and it cannot be neglected. 

As a comment let us notice that knowing the second 
moment of P{v) jlj] enables us to calculate the effec- 
tive temperature T of our systems defined as an averaged 
square velocity 
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To check the validity of our calculations we also made 
Monte Carlo simulations of this model. Results, which 
are shown in Fig. |^ confirm the bilateral exponential de- 
cay of P{v) although the accuracy is this time much 
lower. Our Monte Carlo data are rescaled in such a 
way that their variance matches that obtained using the 
Fourier transform method and a very good agreement is 
seen even on the logarithmic scale. 

In conclusion, we have shown that the Maxwell model 
proposed by Ben-Naim and Krapivsky has velocity dis- 



tribution decaying as a bilateral exponential. Together 
with the recent results by Ernst and Brito ||l^ it indicates 
that such a decay might be of more generic nature for 
this class of models. As a possible extension, it would be 
desirable to examine some other Maxwell models where 
velocities are not scalars but rather d-dimensional vec- 
tors. Actually, such models were already studied and the 
analysis indicates that, for increasing d, correlations be- 
tween velocities and deviation from the pure Gaussian 
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FIG. 2: Velocity distribution P{v) calculated using Monte 
Carlo simulations. Simulations were made for A'^ = 10^ par- 
ticles. Continuous lines are the results obtained from Fourier 
inversion. The inset shows our data in the semi-logarithmic 
scale. 



distribution decrease ||l3|. One possibility is that there 
might be a critical dimension dc and such that for d < dc 
the velocity distribution has a bilateral exponential de- 
cay (as in the present model) while it it has Gaussian 
decay for d> dc. Analysis of such models is however left 
for the future. 
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